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Abstract  Experimental  studies  have  demonstrated  that  botulinum  neurotoxin 
serotype  A  (BoNT/A)  causes  flaccid  paralysis  by  a  multi-step  mechanism.  Follow¬ 
ing  its  binding  to  speciflc  receptors  at  peripheral  cholinergic  nerve  endings,  BoNT/A 
is  internalized  by  receptor-mediated  endocytosis.  Subsequently  its  zinc-dependent 
catalytic  domain  translocates  into  the  neuroplasm  where  it  cleaves  a  vesicle-docking 
protein,  SNAP-25,  to  block  neurally  evoked  cholinergic  neurotransmission.  We  tested 
the  hypothesis  that  mathematical  models  having  a  minimal  number  of  reactions  and 
reactants  can  simulate  published  data  concerning  the  onset  of  paralysis  of  skeletal  mus¬ 
cles  induced  by  BoNT/A  at  the  isolated  rat  neuromuscular  junction  (NMJ)  and  in  other 
systems.  Experimental  data  from  several  laboratories  were  simulated  with  two  differ¬ 
ent  models  that  were  represented  by  sets  of  coupled,  flrst-order  differential  equations. 
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In  this  study,  the  3-step  sequential  model  developed  by  Simpson  (J  Pharmacol  Exp 
Ther  212:16-21,1980)  was  used  to  estimate  upper  limits  of  the  times  during  which 
anti-toxins  and  other  impermeable  inhibitors  of  BoNT/A  can  exert  an  effect.  The 
experimentally  determined  binding  reaction  rate  was  verified  to  be  consistent  with 
published  estimates  for  the  rate  constants  for  BoNT/A  binding  to  and  dissociating 
from  its  receptors.  Because  this  3 -step  model  was  not  designed  to  reproduce  tempo¬ 
ral  changes  in  paralysis  with  different  toxin  concentrations,  a  new  BoNT/A  species 
and  rate  (ks)  were  added  at  the  beginning  of  the  reaction  sequence  to  create  a  4- step 
scheme.  This  unbound  initial  species  is  transformed  at  a  rate  determined  by  to  a  free 
species  that  is  capable  of  binding.  By  systematically  adjusting  the  values  of  ks,  the 
4-step  model  simulated  the  rapid  decline  in  NMJ  function  (ks  >  0.01),  the  less  rapid 
onset  of  paralysis  in  mice  following  i.m.  injections  (ks  =  0.001),  and  the  slow  onset 
of  the  therapeutic  effects  of  BoNT/A  (ks  <  0.001)  in  man.  This  minimal  modeling 
approach  was  not  only  verified  by  simulating  experimental  results,  it  helped  to  quan¬ 
titatively  define  the  time  available  for  an  inhibitor  to  have  some  effect  (finhib)  and  the 
relation  between  this  time  and  the  rate  of  paralysis  onset.  The  4-step  model  predicted 
that  as  the  rate  of  paralysis  becomes  slower,  the  estimated  upper  limits  of  (finhib)  for 
impermeable  inhibitors  become  longer.  More  generally,  this  modeling  approach  may 
be  useful  in  studying  the  kinetics  of  other  toxins  or  viruses  that  invade  host  cells  by 
similar  mechanisms,  e.g.,  receptor-mediated  endocytosis. 
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Abbreviations 

BoNT/A  Botulinum  neuro toxin  type  A 
NMJ  Neuromuscular  junction 

tio  Time-to- 10%  peak  tension 

tinbib  Time  at  which  >  10%  of  unbound  BoNT/A  is  present  and  available 

for  an  inhibitor  to  exert  some  effect 

Introduction 

The  neurotoxins  from  Clostridium  hotulinum  and  several  related  species  represent 
some  of  the  most  lethal  substances  known  [1-3].  The  signs  and  symptoms  include 
flaccid  paralysis  of  the  voluntary  muscles,  respiratory  distress  and  death.  The  onset 
times  and  durations  of  paralysis  depend  on  the  serotype  involved,  the  exposure  route 
and  the  intoxicating  dose.  As  summarized  in  [4],  the  public  is  becoming  increasingly 
aware  of  the  roles  of  botulinum  neurotoxins  as  food  poisoning  agents,  as  potential 
bioweapons  [1,2, 5, 6],  and  as  approved  treatments  for  various  neurologic  indications 
and  other  clinical  uses  [7].  Significant  resources  [8,9]  have  been  devoted  to  the  large- 
scale  production  of  heptavalent  botulism  antitoxin  [10].  Complementary  research  to 
engineer  and  develop  high-affinity,  monoclonal  neutralizing  antibodies  is  also  being 
conducted  [11]. 
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The  bacteria  express  these  toxins  as  single  chain  polypeptides  (MW  ~  150kDa) 
which  are  later  post-translationally  modified  to  form  two  chains  (heavy,  lOOkDa  and 
light,  50kDa)  that  are  covalently  linked  by  a  disulfide  bridge.  The  C-terminal  half  of 
the  heavy  chain  specifically  binds  to  extracellular  acceptors  at  peripheral  cholinergic 
nerve  terminals  [12]  that  innervate  striated  and  smooth  muscles.  A  process  resem¬ 
bling  receptor-mediated  endocytosis  internalizes  the  toxin-bound  receptor.  As  the 
intravesicular  environment  becomes  acidic  (pH  5),  the  N-terminal  half  of  the  heavy 
chain  helps  form  cation-selective  channels  that  may  be  involved  in  allowing  the  escape 
of  the  toxic  moiety  (presumably  the  catalytic  light  chain  or  its  derivatives)  into  the 
neuroplasm  (reviewed  in  [13]).  The  toxic  fragment  is  a  zinc-dependent  protease  that 
cleaves  at  distinct  sites  and  in  a  serotype-specific  manner  one  or  more  of  the  SNARE 
proteins  (SNAP-25,  syntaxin  and  VAMP)  involved  in  the  synaptic  vesicle-mediated 
release  of  acetylcholine.  Once  internalized,  BoNT  is  no  longer  susceptible  to  circulat¬ 
ing  neutralizing  antibodies  or  other  impermeable  inhibitors  of  its  toxicity.  This  homol¬ 
ogous  family  of  proteins  are  grouped  into  seven  immunologically  distinct  serotypes 
(BoNT/A-G)  [3, 14].  SNAP-25  is  cleaved  by  BoNT  serotypes  A,  E  and  Cl,  syntaxin 
is  cleaved  by  BoNT/Cl,  and  VAMP  is  cleaved  by  the  remaining  BoNT  serotypes  [14]. 

The  present  study  was  designed  to  extend  a  data-driven  minimal  model  devel¬ 
oped  by  Simpson  [15]  that  described  the  kinetics  of  botulinum  neurotoxin  serotype  A 
(BoNT/A)  at  the  neuromuscular  junction  (NMJ)  in  producing  paralysis  in  vitro.  This 
original  deterministic  model  consisted  of  a  sequence  of  reactions  based  on  the  known 
mechanism  of  BoNT/A  action,  namely,  binding  to  specific  receptors  located  at  cholin¬ 
ergic  nerve  terminals,  translocating  into  the  neuroplasm  and,  in  turn,  exerting  a  toxic 
effect.  All  three  steps  were  separately  examined  experimentally  and  quantitatively 
characterized  by  apparent  first-order  reaction  rates.  Modifications  were  introduced 
in  our  study  to  allow  for  the  changes  in  paralysis  time  course  seen  under  different 
in  vivo  conditions  [16-18].  We  also  developed  a  quantitative  relationship  between  the 
onset  rate  of  paralysis  and  the  time  that  is  available  to  neutralizing  antitoxins  or  other 
non-permeable  countermeasures  to  exert  some  inhibitory  effect. 


Methods 

Reaction  rates  and  other  rate  constants 

The  experimentally  measured  reaction  rate  of  BoNT/A  binding  to  the  in  vitro  NMJ 
preparation  model  was  compared  to  previously  determined  association  and  dissociation 
(on  and  off)  rate  constants  reported  in  [19]  for  rat  brain  synaptosomes.  These  micro¬ 
scopic  rate  constants  were  adjusted  for  a  ~20°C  difference  in  temperature  because 
the  NMJ  experiments  were  conducted  at  35°C  and  the  synaptosomal  studies  were 
performed  at  4°C.  The  experimentally  determined  value  for  the  temperature  coeffi¬ 
cient,  Qio,  for  the  binding  step  was  used  as  the  multiplicative  factor  by  which  a  rate 
constant  is  increased  when  increasing  the  temperature  by  10°C  [20].  It  was  assumed 
that  the  measured  Qio  value  was  a  constant  for  this  20° C  range  of  temperatures,  that  a 
steady  state  of  equilibrium  was  achieved  and  that  toxin  was  not  internalized  following 
binding. 
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Reaction  equations:  3 -Step  Model 

The  present  models  for  the  receptor  binding,  the  nerve-terminal  entry  and  the 
subsequent  toxic  (lytic)  activity  of  BoNT/A  were  developed  using  the  same  general 
reaction  schemes  and  numerical  procedures  that  were  previously  used  for  modeling 
the  kinetics  of  neurotransmission  at  peripheral  and  central  synapses  [2 1 , 22] .  Reactions 
describing  the  3-step  model  were  derived  from  [15]  for  BoNT/A  intoxication  at  the 
in  vitro  rat  NMJ: 


free  bound  ^  trans  ^  lytic 

where  the  system  of  coupled  first-order  differential  equations  is: 

d\free]ldt  =  —kslfree] 
d[bound]/dt  =  kslfree]  —  kj  [bound] 
d[trans]/dt  =  kr  [bound]  —  kiUrans] 
d[lytic]/dt  =  kiUrans]. 


(I) 


(la) 

(lb) 

(lc) 

(ld) 


This  3 -step  model  represents  a  minimal  version  of  the  original  model  for  a  free  amount 
of  neurotoxin  to  bind,  translocate  into  the  neuroplasm  and  generate  lytic  activity  [15]. 
Each  species  of  BoNT/A  is  associated  with  a  different  environment:  extracellular  in 
solution  (free),  on  the  surface,  bound  to  a  receptor  (bound),  intracellular  endocytotic 
vesicle  (trans),  and  intracellular  neuroplasm  (lytic).  The  three  irreversible  steps  even¬ 
tually  lead  to  the  blockade  of  vesicle-mediated  neurotransmission  and  the  development 
of  complete  muscle  paralysis  in  the  NMJ  experiments  [15] .  The  binding  step  implicitly 
contains  diffusion  of  free  BoNT/A  (0.1  nM)  to  its  binding  sites  at  synaptic  termini. 
The  bound  toxin  undergoes  translocation  into  the  neuroplasm.  Intracellularly,  the  toxic 
moiety  exerts  it  enzymatic  proteolytic  activity  on  SNAP-25  [14].  Embedded  in  this 
species  is  the  paralysis  that  occurs  postsynaptically.  Rates  of  reaction  were  taken  from 
[15]  where^5  =  0.058min“^  kr  =0.141min“^  andkr  =  0.013 min“^  This  minimal 
3-step  scheme  is  a  simplification  from  [15]  because  it  omits  three  unknown  parameters 
representing  binding  sites  for  each  species  of  toxin.  An  analytic  expression  derived 
for  the  lytic  species. 


lytic(t)  =  [T]tot 

[e~^'^^kBkL  (kB  -  ki)  -  e~^^^kBkT  (kB  -  kr)  +  e~^^^kLkT  (ki  -  kr)) 
(kB  —  ki)  (kB  —  kr)  (ki  —  kr) 


(- 


)■ 

(le) 


was  used  to  verify  the  selection  of  the  time  steps  for  the  numerical  solution  and 
to  establish  further  confidence  in  our  minimal  approach.  [T]tot  represents  the  total 
normalized  amount  of  toxin.  Analytic  expressions  for  the  other  species  are  available 
from  the  Electronic  Supplementary  Material. 
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Reaction  equations:  4-Step  Model 

The  4-step  model  is  an  extension  of  the  3-step  model.  BoNT/A  is  present  at  t  =  0  in 
the  bulk  solution  (S)  and  the  rate  constant,  ks,  simulates  the  rate  of  movement  of  the 
non-binding  species  of  BoNT/A  from  the  bulk  medium  to  the  extracellular  solution 
that  surrounds  the  nerve  termini: 


bulk  ^  free  bound  ^  trans  ^  lytic 

(ID 

d[bulk]/dt  =  —ks[bulk] 

(2a) 

d[free]ldt  =  ks[bulk]  —  kB[free]. 

(2b) 

In  the  3-step  model,  initial  conditions  (t  =  0)  for  the  normalized  molecular  species 
of  BoNT  were  [free]  =  1.0,  [bound]  =  0  and  [lytic]  =  0,  while  in  the  4-step  model  these 
species  were  initially  set  to  zero  and  the  [bulk]  =  1 .0.  As  in  the  3-step  model,  an  analytic 
expression  was  derived  for  the  lytic  species: 


lytic(t)  =  [T]tot  - 


[T]tot 


/  e~^'^^kBkLks  (kB  -  ki)  (kB  -  ks)  (ki  -  ks)  -  \ 
e^^^^kBkikj  {kB  —  ki)  {kB  —  kr)  {ki  —  kr)  + 
e~^^^kBkskT  {kB  —  ks)  {kB  —  kr)  {ks  —  kr)  — 
\e-^Btk^kskT  {ki  -  ks)  {ki  -  kr)  {ks  -  kr)  / 


{kB  —  ki)  {kB  —  ks)  {ki  —  ks)  {kB  —  kr)  {ki  —  kr)  {ks  —  kr) 

(2c) 


Analytic  expressions  for  the  other  species  in  this  model  are  available  from  the  Elec¬ 
tronic  Supplementary  Material.  Values  of  ks  were  systematically  changed  from  0.0001 
to  15  min“^ .  In  these  models  of  in  vivo  systems  there  is  no  explicit  clearance  (i.e.,  loss 
of  toxin  from  those  systems),  only  an  implicit  one  that  is  embedded  in  the  k^  rate  value. 
It  was  assumed  that  the  amount  of  BoNT  that  was  injected  intramuscularly  became 
distributed  in  the  extra- vascular  space,  and  was  partially  taken  up  by  the  circulatory 
system  or  in  various  organs  (e.g.,  liver,  spleen,  kidney)  as  has  been  demonstrated  in 
[4]  following  i.v.  injections  of  BoNT.  The  kinetic  parameters  used  in  these  minimal 
models  are  tabulated  at  http://botdb.abcc.ncifcrf.gov. 

Diaphragm  muscle  function  (peak  twitch  amplitude  induced  by  electrical  stimula¬ 
tion  of  the  phrenic  nerve)  was  defined  as  the  normalized  amount  of  tension  remaining 
after  the  formation  of  the  lytic  species: 

tension  =  1.0  —  [lytic].  (3) 

Values  of  tio  are  the  times  at  which  the  model  predicts  that  the  peak  twitch  tension 
is  10%  of  the  maximum  tension,  i.e.,  a  90%  reduction  in  tension.  This  level  of  residual 
tension  was  used  in  these  in  vitro  experiments  because  complete  paralysis  is  difficult  to 
accurately  measure  [15].  The  twitch  responses  were  measured  after  exposure  to  300  ng 
( 10“  M)  BoNT/A  in  20  ml  bathing  chambers  as  in  previous  in  vitro  NMJ  studies  [23] . 
Nerve-muscle  preparations  were  exposed  to  toxin-containing  media  that  were  made 
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with  0.5  mM  and  5.0 mM  to  prevent  toxin  internalization.  After  a  60 min 

incubation,  the  preparations  were  washed  with  saline  then  soaked  in  control  medium. 

Other  response  measures  were  used  in  the  mouse  in  vivo  studies  [7,17]  and  in 
clinical  studies  [18].  In  [7],  time  spent  by  mice  on  an  exercise  wheel  was  monitored 
during  the  onset  and  recovery  of  the  BoNT/A-induced  paralysis  (0.2-1.5ng).  The 
study  in  [17]  measured  toe  abduction  reflexes  as  the  functional  parameter  affected  by 
the  i.m.  injection  of  ~26ng  of  BoNT/A.  In  the  clinical  study  [18],  the  patient  reported 
the  alleviation  of  subjective  symptoms  of  cervical  dystonia  following  i.m.  injections  of 
BoNT/A  which  are  1-4  ng  according  to  [24]  for  this  neurologic  disorder.  In  all  of  the 
simulations  using  these  in  vivo  data,  the  term  “response”  replaces  “tension”  in  Eq.  3. 

In  this  study,  we  estimated  the  upper  limit  of  time  that  non-permeable  inhibitors  of 
BoNT/A  (e.g.,  an  antibody)  have  some  effect.  The  value  of  this  available  time  (tmhib) 
was  operationally  defined  as  the  time  at  which  the  amount  of  unbound  BoNT/A  is  at 
least  10%  of  the  initial  unbound  species.  We  are  using  this  metric  from  the  point  of 
view  that  it  is  the  toxin  that  needs  to  be  inhibited  by  an  antibody  or  other  non-permeable 
inhibitor.  Because  finhib  is  quantitative,  it  can  be  later  refined  when  the  time  course  of 
the  unbound  amount  of  BoNT/A  is  experimentally  resolved.  In  the  3 -step  model,  the 
relative  value  of  the  unbound  BoNT/A  is  equal  to  the  value  of  the  free  species,  while 
in  the  4-step  model,  the  relative  value  of  the  unbound  BoNT/A  is  equal  to  the  sum  of 
the  values  of  the  free  and  the  bulk  species. 

Computer  programs  and  statistical  methods 

The  kinetic  schemes  were  visually  constructed  and  edited  interactively  using 
JDesigner  which  is  part  of  the  open  source  Systems  Biology  Workbench  (SBW)  [25]. 
This  program  saves  models  in  SB  ML  format  so  that  they  can  be  used  as  a  front-end 
software  component  to  any  SBW-compatible  equation  solver.  The  ordinary  first-order 
differential  equations  describing  the  reactions  built  with  JDesigner  were  solved  by 
Jarnac,  a  script  language-based  program  for  analyzing  biochemical  reaction  networks 
[25].  Mathematica(g)  (Wolfram  Research  Inc,  Champaign,  Illinois)  generated  analytic 
expressions  that  were  used  to  verify  the  selection  of  the  time  steps  for  the  numeri¬ 
cal  solutions  by  Jarnac  and  to  establish  further  confidence  in  our  minimal  approach. 
A  spreadsheet-compatible  data-analysis  utility  was  used  to  digitize  and  store  two- 
dimensional  coordinates  of  data  points  from  selected  figures  shown  in  [4,7,15-17] 
and  [18]  (Datatrend  Software,  Raleigh,  North  Carolina).  Linear  regression  analy¬ 
sis  of  the  form  y  =  yo  +  ax  was  performed  on  the  data  in  Lig.2  to  calculate  coeffi¬ 
cients  of  determination  (R^)  and  values  of  the  slope  (a)  and  intercept  (yo)  along 
with  standard  errors  of  the  mean  (S.E.M.)  using  SigmaPlot  (ver.  9.01,  2004,  Sys- 
tat  Software,  Inc.,  Chicago,  Illinois).  Area  under  the  curve  (AUC)  calculations  were 
performed  on  the  numerical  outputs  (Lig.  1)  by  summing  the  amount  of  free  toxin 
(ng  units)  present  at  each  time  increment  (60 min)  in  the  3-step  model.  Because 
of  the  small  amount  of  free  toxin  predicted  for  the  simulation  in  Lig.  4,  only  the 
amount  of  bulk  toxin  was  summed  at  each  time  increment  (1  min)  in  the  4- step  model. 
Elimination  rates  of  the  unbound  (i.e.,  free  or  bulk)  toxin  species  were  calculated 
by  taking  the  amount  (ng  units)  at  t  =  0  and  dividing  by  the  AUC  value  (min-ng 
units). 
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Results 

Relating  the  binding  reaction  rate  to  the  receptor  binding  and  dissociation 
rate  constants  for  BoNT/A 

Before  examining  the  3 -step  model,  we  tested  the  hypothesis  that  the  receptors  for 
BoNT/A  are  not  saturated  at  0.1  nM,  the  concentration  used  in  [15].  In  addition,  we 
evaluated  the  relation  between  the  measured  apparent  first-order  binding  reaction  at 
the  rat  NMJ  (^5)  [15]  and  the  association  ikon)  and  dissociation  {koff)  microscopic 
rate  constants  that  have  been  measured  using  rat  cerebrocortical-derived,  synaptosome 
preparations  [19].  The  values  for  kon  =  0.0138 nM“^min“^  and  koff  =  0.0072min“^ 
at  4°C  [19]  correspond  to  sl  Kd  =  0.52nM.  [15]  the  rate  of  paralysis  in  the  rat  NMJ 
preparation  (ks  =  0.058  min“^  was  measured  under  conditions  of  low  [Ca^+],  high 
[Mg^+]  containing  solutions  at  35-36°C  to  block  the  internalization  step  (Fig.  2  of 
[  1 5 ] )) .  In  keeping  with  experimental  observations  [  1 5 ] ,  it  was  assumed  that  low  [Ca^+ ] , 
high  [Mg^+] -containing  solutions  did  not  affect  these  on  and  off  rate  constants  in  this 
reversible  reaction: 


U  +  R  ^  B  (III) 

^Off 

where  [U]  =  [BoNT/A] (unbound),  [R]  =  binding  site  (productive  ecto-acceptor) 
[26],  and  [B]  =  [BoNT/A] (bound).  Using  the  experimentally  determined  value  for 
Qio  =  1.6  for  the  binding  step  [15],  the  temperature-adjusted  values  for  kon  and  koff 
at  35°C  were  estimated  to  be  0.056  nM“^  min“^  and  0.029  min“^  respectively,  and 
are  used  in  the  equations  below.  The  total  toxin  and  binding  site  concentrations  at 
t  >  0  are 


[T]tot  =  [U]  +  [B];[R]tot  =  [R]  +  [Bl 


(4) 


The  binding-site  occupancy  was  assumed  to  be  described  by  the  mass-law  formu¬ 
lation  to  give  an  equilibrium  dissociation  constant 


koff  ^  [U][R] 
kon  [^] 


(5) 


in  which  the  species  concentrations  are  substituted  for  their  activities.  These  terms 
are  related  to  the  experimentally  measured  value  for  an  apparent  first-order  rate 
constant,  that  is  the  reciprocal  of  the  relaxation  time  constant  (r)  for  a  reversible, 
single-step  bimolecular  reaction  [27]: 


ks  —  l/”^  —  i[T]tot  ”1”  [^]tot)kon  H”  ^off  (6) 

Combining  and  rearranging  Eqs.  4-6  yields  a  quadratic  expression 

[/?]2  +  {2Kd  -  kB/kon}[R]  +  (KAT],o,  -  Kdke/kon  +  Kf=0  (7) 
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in  which  the  positive  root  of  [R]  is: 


[R] 


-{IKd  -  kB/Kn]+^{lKd  -  kB/Kn)}'^+A{Kdmtot 

2 


KdkB/kon  +  Kj) 


(8) 


Solving  for  [R]  allows  values  for  [U]  and  [B]  to  be  calculated 

^  ([T]tot/Kd) 

[7?]  +  (1/Kd) 

[B]  =  [T]tot  -  [Ul 


(9) 

(10) 


Using  this  approach  with  a  dimensionless  value  of  Ttot  =  1-0,  the  values  for  the 
reaction  species  at  35°C  (assuming  a  steady  state  of  equilibrium  and  no  toxin  inter¬ 
nalization)  were  calculated  to  be  [R]  =  8.6,  [U]  =0.69  and  [5]  =  0.31.  From  these 
results,  the  values  of  and  the  corresponding  on  and  off  rate  constants  support  the 
idea  that  the  binding  sites  for  BoNT/A  are  not  saturated.  These  results  also  support  the 
hypothesis  that  the  binding  reaction  rate  for  BoNT/A  measured  in  the  isolated  NMJ 
at  35-36°C  can  be  feasibly  compared  with  the  microscopic  on/off  rate  constants  for 
BoNT/A  measured  in  rat  brain  at  4°C. 


The  3 -step  model:  estimates  of  tinhib  and  the  frequency  dependence  of  toxicity 
at  the  NMJ 

As  depicted  in  Fig.  1,  the  experimental  data  collected  from  the  in  vitro  NMJ  (adapted 
from  Fig.  4  of  [15])  was  simulated  with  a  minimal  3 -step  model.  In  those  experiments, 
BoNT/A  was  incubated  in  low-Ca^^ ,  high-Mg^+  containing-media  for  40  min  to  allow 
near-maximal  binding.  A  second  incubation  period  with  normal  media  allowed  for 
near-maximal  translocation  ("^90%  in  16  min)  which  was  followed  by  the  commence¬ 
ment  of  phrenic  nerve  stimulation  to  monitor  the  gradual  decline  in  twitch  tension. 
In  contrast,  BoNT/A  was  applied  in  this  and  the  4-step  model  at  t  =  0.  The  time 
courses  of  the  four  BoNT/A  species  can  be  visually  compared  during  the  development 
of  BoNT/A-induced  paralysis  (time-to-10%  tension,  tio  =  204min).  The  numerical 
solution  generated  by  Jarnac  for  the  tension  curve  using  a  time  step  of  1  min  (Fig.  1, 
black  line)  was  verified  with  an  analytic  expression  (white  curve  superimposed  on 
black  line,  Eqs.  le  and  3).  The  time  course  of  the  free  BoNT/A  species  provides  a 
quantitative  metric  for  estimating  the  time  available  for  an  inhibitor  to  exert  some 
effect  (tinhib),  as  defined  in  the  Methods.  In  this  3-step  model,  tio  was  predicted  to  be 
five  times  longer  than  finhib-  This  model  predicted  that  to  exert  some  degree  of  inhi¬ 
bition,  impermeable  neutralizing  antibodies  or  BoNT/A  inhibitors  need  to  be  applied 
within  40  min  of  the  addition  of  0.1  nM  BoNT/A  to  the  bathing  chamber,  i.e.,  before 
a  ~20%  reduction  of  nerve-evoked  twitch  tension  occurs. 

The  use-dependence  associated  with  BoNT/A  at  the  NMJ  (data  adapted  from  Fig.  5 
of  [15])  is  illustrated  in  Fig.  2A  in  which  paralysis  times  decreased  (smaller  values  of 
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Fig.  1  Simulation  of  experimental  NMJ  data  with  the  3-step  model.  The  time  courses  of  the  model  species 
can  be  compared  with  the  development  of  BoNT/A-induced  paralysis  (filled  symbols,  data  from  [15]).  The 
time- to- 10%  tension,  tip  ,  was  204  min.  The  upper  limit  for  the  amount  of  time  available  for  an  inhibitor 
of  BoNT/A  to  exert  some  effect,  finhib^  was  <40 min  (see  Methods).  All  BoNT/A  species  are  represented 
as  relative  values.  Dotted  line:  free  BoNT/A;  long  dashes:  bound  toxin;  short  dashes:  translocated  species; 
black  line:  amplitudes  of  peak  twitch  tension  numerically  calculated  from  the  3-step  model;  white  line 
(superimposed  on  black  line):  amplitudes  of  peak  twitch  tension  calculated  from  Eq.  le 


tio)  as  the  frequency  of  stimulating  the  phrenic  nerve  increased.  Results  of  the  linear 
regression  analysis  (see  Methods)  were: 

?10  =  94.6(5.4,  S.E.M  )  -  85.4(35.9,  S.E.M.  )frequency  (R^  =  0.73). 


Changing  the  values  of  either  ks  ox  kj  could  only  partially  simulate  the  experimental 
data.  In  contrast,  altering  the  values  for  the  third  reaction  rate  {ki)  could,  for  the  first 
time  to  our  knowledge,  simulate  the  entire  set  of  data  associated  with  the  frequency- 
dependent  enhancement  of  BoNT/A-induced  paralysis.  Values  for  tio  (Fig.  2B)  overlap 
with  the  experimental  data  points  that  reside  within  a  linear  portion  of  the  predicted 
non-linear  relation  between  values  of  ki  and  normalized  times-to-10%  tension.  The 
values  of  ki  (0.0 13-0.01 9  min“^)  correspond  to  values  of  stimulation  frequencies  that 
are  within  practical  experimental  limits  (Fig.  2C): 

frequency  =  0.54(0.67,  S.E.M. )  +  43.25(4.16,  S.E.M.  )kL  {R^  =  0.98). 


The  4-step  model:  in  vitro  NMJ  versus  in  vivo  systems  and  estimates  for  finhib 

Altering  the  three  original  reaction  rates  {kB^kr  and  k^)  failed  to  produce  the  slower 
onset  rates  of  BoNT/A-induced  paralysis  observed  in  vivo.  Thus,  the  3-step  model 
examined  above  was  extended  and  used  to  simulate  BoNT/A-induced  responses  of 
patients  when  administered  as  an  approved  therapeutic  product  [18]  and  of  mice  when 
injected  as  a  local  paralysant  [7, 17].  In  this  model,  a  new,  initial  reaction  rate  {ks)  was 
created  that  represents  the  movement  of  BoNT/A  (bulk  species)  from  distal  locations 
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Fig.  2  Different  values  of  in 
the  3 -step  model  simulated  the 
frequency-dependence  of  the 
rate  of  onset  of  BoNT/A-induced 
paralytic  activity.  Panel  A:  plot 
of  data  points  where  values  for 
1 10  are  expressed  as  normalized 
values.  Panel  B:  The  3-step 
model  (dashed  line)  produced  a 
non-linear  relation  between 
k/^(0.01  —  0.3min“^)  and  the 
time- to- 10%  tension,  tio-  Within 
the  range  of  experimental  data, 
the  plot  was  linear  and  at  higher 
values  of  approached  an 
asymptotic  value  for  the 
normalized  values  of  tio-  Panel 
C:  Values  of  ki  were  adjusted  to 
simulate  normalized  tio  values 
at  each  of  the  four 
experimentally  used  nerve 
stimulation  frequencies.  Circles: 
experimental  data  from  Fig.  4  in 
[15];  lines  in  Panels  A  and  C: 
least-squares  fits  drawn  through 
the  data  points 


kL(min  *) 


from  its  receptors  to  the  proximal  volume  where  the  BoNT/A  (free  species)  is  available 
for  binding. 

Slow  onset  kinetics  for  the  beneficial  effects  of  treatments  using  BoNT/A  have 
been  published  in  the  clinical  literature.  Normalized  clinical  data  represented  in 
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Fig.  3  Different  values  of  ks  in  the  4-step  model  simulated  the  onset  of  BoNT/A-induced  effects  in  a 
variety  of  in  vitro  and  in  vivo  data  sets.  The  onset  of  paralysis  with  bath-applied  BoNT/A  at  the  rat  NMJ 
(circles,  from  Fig.  1)  is  rapid  and  is  associated  with  ks  values  >0.01  min“^.  In  vivo  data  (time  using  an 
exercise  wheel,  see  Methods)  from  locally  injected  mice  (stars,  data  from  [17];  upward  triangles,  data 
from  [7])  were  simulated  using  ks  =  0.001  min“^ .  Clinical  data  (subjective  measure  of  relief  of  neurologic 
symptoms,  see  Methods;  downward  triangles)  are  from  a  patient  whose  first  session  of  BoNT/A  treatment 
was  on  day  0  (data  from  [18]).  The  maximum  BoNT/A-induced  effect,  achieved  by  day  11,  was  simulated 
with  ks  =  0.00015  min“^ 


Fig.  3  are  from  a  patient  with  cervical-dystonia-related  complaints  receiving  BoNT/A 
intramuscularly  (data  from  Fig.  1,  injection  series  1  in  [18]).  Starting  at  approximately 
day  11,  the  maximum  relief  from  symptoms  was  noted  by  the  patient  (tio=  10.6  days). 
The  rate  constant  /c^was  systematically  varied  from  0.0001  to  10min“^  to  provide  an 
estimated  value  of  0.00015  min“^  for  this  model. 

Other  examples  from  two  separate  experimental  studies  using  mice  are  also  included 
in  Fig.  3  [7,17].  Both  sets  of  data  following  local  intramuscular  injections  could  be 
simulated  (ks  =  0.001  min“^).  These  data  show  a  more  rapid  onset  than  the  clinical 
data  but  are  slower  than  the  onset  of  paralysis  seen  with  the  in  vitro  rat  NMJ  data  (see 
also  Fig.  1).  In  this  4-step  model,  the  NMJ  data  could  be  simulated  with  values  of  ks  > 
O.Olmin-i. 

A  comparison  was  also  made  of  the  time  courses  for  the  different  toxin  species  in 
the  3-step  model  for  the  NMJ  and  the  4-step  model  for  the  in  vivo  data  examined  here. 
As  for  the  model  tension  prediction  in  Fig.  1,  the  numerical  solution  for  the  tension 
curve  (Fig.  4,  black  line)  was  verified  analytically  (Fig.  4;  white  curve,  Eqs.  2c  and  3). 
As  mentioned  above,  the  3-step  NMJ  model  predicted  that  the  time-to-10%  tension 
(tio)  was  5  times  slower  than  finhib-  In  contrast,  the  clinical  data  had  a  much  slower 
onset  (tio  =  15,  450min  or  10.6 days)  for  the  BoNT-induced  response,  whose  value 
was  almost  identical  to  that  of  finhib  (Fig-  4).  In  this  case,  unbound  BoNT/A  was  com¬ 
posed  largely  of  the  bulk  species  with  a  negligible  presence  of  the  free  BoNT/A  species. 
Thus,  in  the  4- step  model  for  the  clinical  data,  the  time  course  of  the  unbound  BoNT/A 
species  was  predicted  to  develop  at  essentially  the  same  rate  as  the  time  course  of 
paralysis. 
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Fig.  4  Time-dependence  of  the  BoNT/A  species  from  the  4-step  model  that  simulated  the  clinical  data. 
The  time  courses  of  the  different  species  can  be  compared  with  the  relatively  slow  loss  of  tension 
(tio  =  15,450min  or  106  days)  of  BoNT/A-induced  paralysis.  In  contrast  to  the  in  vitro  NMJ  (Fig.  1),  the 
time  courses  of  the  unbound  BoNT/A  and  the  lytic  species  overlapped,  making  the  estimated  upper  limit 
value  of  iinhib  longer  in  vivo  than  at  the  isolated  NMJ.  Closed  circles:  subjective  clinical  data  from  [18]; 
black  line:  responses  were  calculated  numerically  from  the  4-step  model  using  ky  =  0.00015  min“^ ;  white 
line:  responses  were  calculated  from  Eq.  2c 


The  area  under  the  curve  (AUC)  calculations  were  done  to  compare  the  temporal 
profiles  of  the  unbound  species  of  BoNT  using  in  vitro  NMJ  and  in  vivo  human  data. 
The  areas  for  unbound  BoNT/A  are  about  the  same  (<20%  difference)  for  the  isolated 
NMJ  and  for  the  in  vivo  injections  in  man.  The  NMJ  was  exposed  to  300  ng  while 
the  patient  received  an  i.m.  injection  of  approximately  1  ng  BoNT/A.  The  amount  of 
unbound  BoNT  at  the  NMJ  decreased  by  95%  in  50 min.  In  man,  the  same  decrease 
occurred  in  about  20,000  min  (400  times  longer).  The  AUC  for  the  NMJ  is  5 174  min-ng 
using  the  4-step  model  with  >  0.01  min“^,  while  the  AUC  for  man  is  6406  min- 
ng  using  the  4-step  model  with  =  0.00015  min“^.  As  the  value  of  k^  decreased, 
the  amount  of  the  unbound  BoNT  at  any  time  point  increased.  Given  the  presumed 
differences  in  toxin  purities  and  the  uncertainty  in  the  amount  the  patient  received, 
these  area  values  are  remarkably  similar.  The  rate  of  disappearance  of  the  unbound 
BoNT/A  is  400  times  faster  at  the  NMJ  but  this  preparation  also  received  about  300 
times  more  BoNT/A  than  the  patient.  Thus,  the  calculated  areas  are  similar  due  to 
offsetting  differences  in  these  rates  and  the  amounts  of  toxin  used. 

To  demonstrate  the  relation  between  the  upper-limit  estimate  for  tinhib  and  time- 
to-10%  tension  (tio),  the  values  of  ks  were  varied  (Fig.  5).  The  symbols  in  this  figure 
are  the  same  as  in  Fig.  3  for  representing  ks  values  from  the  NMJ,  from  the  i.m. 
injections  in  mice  and  from  the  clinical  data.  The  4- step  model  predicted  that  as  ks 
decreased,  the  absolute  values  for  tio  and  finhib  increased.  The  relations  between  ks  and 
the  upper-limit  estimates  of  finhib  and  times-to-10%  tension  were  non-linear  (Fig.  5 A). 
The  largest  component  of  the  unbound  BoNT/A  species  was  predicted  to  be  the  bulk 
species  for  the  data  from  the  mouse  local  injections  and  from  the  clinical  data  where 
ks  <  0.001  min“^  These  calibration  curves  also  provide  upper-limit  estimates  for 
tinhib  that  Can  be  directly  related  to  the  duration  for  which  a  pre-defined  threshold  of 
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Fig.  5  Effects  of  varying  ky  on  the  upper-limit  estimates  of  hnhib  onset  times  for  BoNT/A-induced 
paralysis.  Panel  A:  the  two  durations,  the  times-to-10%  peak  tension  (tio)  and  the  hnhib^  increase  and  con¬ 
verge  as  the  values  of  ks  decrease.  Panel  B:  the  predicted  concentration  of  unbound  BoNT/A  (decomposed 
into  bulk  and  free  species)  at  hnbib^  is  dependent  on  the  value  of  ky.  The  placement  of  symbols  along  the 
x-axis  in  both  panels  (upward  triangles,  downward  triangles,  stars,  and  circles)  are  defined  in  Fig.  3  and 
depict  which  of  the  unbound  species  and  values  of  ks  are  associated  with  those  simulations.  For  clarity, 
these  symbols  were  placed  just  above  the  x-axis  to  correspond  with  the  appropriate  values  of  k^- 

tension  is  sustained.  The  plots  in  Fig.  5B  over  a  range  of  values  of  ks  show  the  relative, 
decomposed  values  for  unbound  BoNT/A  at  tmhib  for  the  bulk  and  free  species.  The 
same  trends  were  obtained  with  this  model  at  other  levels  of  post-exposure  tension 
(50  and  90%)  that  were  examined  (data  not  shown). 


Discussion 

The  results  of  this  study  support  the  hypothesis  that  a  minimal  approach  can  model 
responses  from  in  vitro  and  in  vivo  systems.  This  common  tactic  of  using  a  reduction¬ 
ist  approach  was  applied  because  computing  the  BoNT/A-induced  perturbations  in 
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nerve-evoked  muscle  twitches  is  beyond  the  practicality  of  developing  fully  detailed 
models  of  the  mechanisms  of  neurotoxin  action  upon  vesicle-mediated  chemical  neu¬ 
rotransmission.  Several  implicit  assumptions  were  made  in  the  two  models  examined. 
Most  importantly,  it  was  assumed  for  computational  convenience  that  the  models  were 
independent  of  compartments  and  that  uniform  distributions  existed  for  all  the  BoNT/A 
species  and  binding  sites.  Despite  these  simplifications  and  our  incomplete  knowl¬ 
edge  of  cellular  structures,  relevant  results  emerged  using  biochemical  observables 
as  guides,  specifically,  the  binding  reaction  rate  [15]  together  with  the  microscopic 
binding  and  dissociation  rate  constants  [19].  These  values  were  consistent  with  classic 
methods  in  ligand-binding  studies  in  which  the  inequality  [R]tot  >  [7^1  tot  is  assumed 
to  exist  within  a  well- stirred  system  that  is  far  from  a  state  of  binding- site  saturation 
[28].  Our  calculations  also  show  that  various  types  of  binding  data  from  different 
tissues  (central  and  peripheral  neurons)  and  temperatures  (4  to  35-36°)  are  related. 

More  insights  were  gained  using  the  minimal  models  about  the  temporal  relations 
between  the  extent  of  paralysis  and  the  time  that  is  available  for  an  inhibitor  to  exert 
some  effect,  finhib,  under  different  in  vitro  and  in  vivo  conditions.  These  simulations 
predicted  upper-limit  estimates  for  finhib  during  which  10%  of  the  maximum  amount 
of  the  unbound  BoNT/A  species  is  available  for  interacting  with  impermeable  antag¬ 
onists.  Perhaps  the  most  important  concept  to  emerge  from  these  minimal  models  is 
that  as  the  onset  rate  of  paralysis  becomes  slower,  i.e.,  as  the  value  of  ks  decreases,  it 
also  becomes  a  better  indicator  of  the  upper-limit  for  when  finhib  and  the  times-to-10% 
tension  (tio)  converge  (Fig.  5).  The  value  of  finhib  is  not  a  static  amount  of  time  but, 
rather,  depends  on  the  rate  of  paralysis  onset. 

These  modeling  results  are  related  to  a  recent  study  of  iatrogenic  botulism.  The 
existence  of  detectable  serum  levels  of  BoNT/A  (12-24  mouse  LD50)  in  three  out  of 
four  case-patients  2-3  days  following  extremely  high  i.m.  doses  (up  to  eight  million 
mouse  LD50)  of  research-grade  toxin  [29]  suggests  that  the  anti-toxin  administration 
may  have  had  a  beneficial  effect.  These  results  also  suggest  that  further  development 
of  animal  models  need  to  be  conducted  in  vivo  to  test  hypotheses  that  use  onset  rates 
of  paralysis  to  obtain  accurate  estimates  of  finhib  for  neutralizing  antibodies  or  imper¬ 
meable,  low  molecular  weight  inhibitors. 

The  binding  step  for  BoNT/A  was  intentionally  simplified  because  of  the  inherent 
complexity  of  the  intoxication  process  and  of  the  variety  of  reactions  comprising 
neurotransmission.  At  least  two  molecular  components  form  the  functional  binding 
site  for  BoNT/A,  a  protein  and  a  ganglioside  [14].  The  ganglioside  is  GTlb  or  a 
related  derivative  [30,31].  Although  the  protein  component  (SV2,  synaptic  vesicle 
protein-2)  has  twelve  putative  transmembrane  spanning  regions  [32,33],  the  location 
of  the  BoNT/A  binding  site  is  on  a  region  of  SV2  that  is  within  a  synaptic  vesicle. 
This  location  implies  that  the  interior  solution  of  this  structure  must  be  in  transient 
contact  with  the  extracellular  medium  so  that  BoNT/A  can  bind  immediately  or  bind 
following  vesicle-pore  closure  and  vesicle  fission  from  the  presynaptic  membrane. 
During  quiescent  periods,  the  majority  of  the  BoNT/A  SV2  binding  sites  are  not  in 
contact  with  the  extracellular  medium  and  are  unavailable  for  binding. 

In  contrast,  with  neural  stimulation,  these  sites  become  exposed  and  available 
for  binding,  and  as  the  frequency  of  stimulation  increases,  the  time-to-10%  tension 
(tio)  decreases  [15]  (Fig.  2).  Although  alterations  in  the  other  two  reaction  rates,  by 
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themselves,  were  insufficient  to  simulate  the  examined  data,  changes  in  all  three  of 
these  reaction  rates  may  be  required  to  accurately  model  the  experimental  data.  Thus, 
the  enhanced  rate  of  paralysis  onset  with  increased  nerve  stimulation  frequency  may 
be  a  complicated  function  of  enhanced  calcium  entry  and  an  enhanced  lytic  effect, 
i.e.,  the  proteolysis  of  the  BoNT/A  substrate  SNAP25. 

Other  simplifications  of  the  binding  process  in  the  4-step  model  involve  the  initial 
step  associated  with  the  reaction  rate,  ks.  As  noted  in  115],  there  are  anatomical 
and  perhaps  other  barriers  that  interfere  with  BoNT/A  readily  diffusing  to  BoNT/A 
binding  sites  in  the  isolated  NMJ.  For  the  simulations  of  both  the  short-term  in  vitro 
and  longer-term  in  vivo  systems,  delays  in  diffusion  to  binding  sites  caused  by  non¬ 
specific  binding  126,34]  also  need  to  be  considered.  Potential  non-specific  binders  of 
BoNT/A  include  serum  proteins  (e.g.,  albumin)  14]  and  extracellularly  located  sites  in 
the  circulatory  system. 

Regarding  the  final  two  steps  in  the  minimal  models,  simplifications  were  made  for 
the  translocation  and  lytic  (toxic)  processes  113, 15,35,36].  The  translocation  step  can 
be  decomposed  into  at  least  a  dozen  steps,  the  most  prominent  being  vesicle  fission 
from  the  plasma  membrane  into  the  neuroplasm,  the  activation  of  energy-dependent 
proton  pumps  with  the  subsequent  decrease  in  intra- vesicular  pH,  the  resultant  con¬ 
formational  changes  to  BoNT/A,  the  reduction  of  its  inter-chain  disulfide  bridge,  the 
ion  channel  formation  by  the  toxin’s  heavy  chain,  and  the  escape  of  the  toxic  moiety 
(presumably  the  50kD  light  chain)  into  the  neuroplasm  where  it  exerts  its  toxic  effect. 
A  simplified  single-step  reaction  is  justified  because  the  rates  associated  with  many 
of  these  reactions  are,  as  yet,  undetermined. 

The  simplified  lytic  process  includes  among  its  most  important  reactions  the  zinc- 
dependent  enzymatic  activity  of  BoNT/A  115]  that  can  be  described  by  Michaelis- 
Menten  kinetics  to  cleave  SNAP-25  137].  This  proteolysis  prevents  the  formation  of 
a  functional  SNARE  protein  complex  (including  synaptobrevin,  syntaxin  and  other 
proteins)  to  block  the  release  of  acetylcholine  and  eventually  block  stimulus-evoked 
muscle  contractions.  The  lytic  step  also  contains  a  steep  non-linear  relation  between  the 
amount  of  substrate  cleaved  and  the  decrease  of  the  evoked  muscle  twitch  amplitude 
138^0]. 

Although  we  selected  a  clinical  study  that  had  an  abundance  of  relevant  quantitative 
data  118],  other  clinical  cases  involving  onset  time  scales  of  days  should  also  be 
computationally  examined  and  compared.  For  example,  regarding  cases  of  food-borne 
14 1  ] ,  iatrogenic  and  wound  botulism  142]  there  is  a  lack  of  extensive  clinical  data  sets  in 
the  publicly  available  literature  and  other  open  sources  to  archive  paralysis  time  courses 
for  future  studies  evaluating  the  temporal  dependence  of  treatment  effectiveness. 

Notwithstanding  the  limitations  of  our  minimal  models,  they  have  allowed  us  to  sim¬ 
ulate  the  onset  kinetics  of  BoNT/A-induced  effects  observed  in  experimental  and  clin¬ 
ical  studies.  Beyond  gaining  a  better  appreciation  of  the  temporal  events  that  underlie 
the  harmful  as  well  as  potentially  beneficial  effects  of  this  family  of  neurotoxins,  these 
models  also  provide  us  with  an  initial,  quantitative  framework  for  making  extrapola¬ 
tions  from  animal  model  to  human  data. 
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